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ABSTRACT 



An initialization of a spectral formulation of the primi 
t i ve equations using a diagnostic divergence is tested for a 
global model. The initial conditions are generated from a 
developing ba roc I i n I ca I I y unstable wave. A sem i - i mp I i c i t ti 
scheme is developed and tested along with the usual exp I icit 
method during the course of the experiments. Results show a 
relatively small effect of a divergent initialization on the 
ensuing integrations. The sem i - i mp I i c i t method shows a ten- 
dency to smooth out high frequency oscillations in local ten 
dene i es . 
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I NTRODUCTION 



The emphasis of this thesis is two-fold. The first point 
is to ga i n insight into the characteristics of a spectral 
primitive equation model. The second is to attempt to take 
advantage of the spectral expansions in an initialization 
techn i que . 

A spectral model uses a Galerkin formulation with the 
spherical harmonics as its choice of basis functions. The 
basis functions are eigenfunctions of the Laplacian operator 
in spherical geometry and, thus, are used to advantage for 
global integrations. An immediate consequence of the Galerkin 
formulation is that spatial derivatives of a single wave are 
computed without the usual truncation error present in finite 
difference methods. This means that the spectral model 
contains no linear phase speed errors. 

The non-linear terms are computed via a transform grid 
in physical space (Orszag, 1971). This procedure is respon- 
sible for making the efficiency of the spectral model com- 
parable with finite difference methods. In addition, the 
non-linear terms may be computed without aliasing. This pre- 
vents non-linear instability and conserves invariant integral 
properties. Physical processes are also computed on the 
transf orm grid. 

The second objective of this thesis is to test a divergent 
initialization for the spectral model. Most static and 
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variational initialization techniques of global primitive 
equation models use a constraint to balance the rotational 
part of the wind with the mass field while setting the 
divergent part to zero. Phillips (I960) has shown that a 
completely non-d i vergent initialization cannot eliminate the 
gravity waves which would be generated and his work suggests 
that the quas i -geostroph i c divergence could tend to suppress 
these gravity wave modes for mid-latitude synoptic scales. 

In addition, if the q ua s i -geo s t ro p h i c divergence had internal 
modes close to those which would be developed by the model, 
such things as the large scale precipitation would be more 
realistic during the first few hours of the model forecast. 

It is an aim of this thesis to develop a diagnostic 
velocity divergence to insert into the fields of an otherwise 
completely balanced system and to test the subsequent effect 
on the spurious gravity wave noise and the developing diver- 
gence. 

A spectral formulation of the primitive equations can 
lend itself to techniques which initialize according to scale 
For example, the q uas i -geostrop h i c divergence might be com- 
puted only for the synoptic scales of the model. Furthermore 
the model conveniently supplies diagnostic information about 
the various scales. 

The final goal of this thesis is to develop a semi-impli- 
cit time differencing scheme and test its effects during the 
first few model hours after initialization. 
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I I . MODEL DESCRIPTION 



The model described here has been developed by Rosmond 
(1977) and the details are presented here for reference 
purposes. Similar formulations have been done by Hoskins 
and Simmons (1975) and by Bourke (1974). The equations for 
an inviscid, adiabatic and hydrostatic atmosphere may be 
wr i tten as 

1.1) 15 = -V*(C+f)V - k • Vx(RTVq + a 1^) 

oT o0 

->-2 

9D ^ • 3V 2 V 

1.2) 15 = f<*Vx(C+f)V - V-(RTVq + a |^) - V^(4> + j ) 

1.3, |9 4 . 79 . „|9 

1.4, |a=.o.v.v, -|2 
1.5) a|i=-RT 

where 

? - vorticity (^ = k • V x ^) 

D - velocity divergence (D = V • V) 

T - temperature 

9 - potential temperature 

IT - surface pressure 

V - horizontal velocity vector 

4> - geopotential height 

R - gas constant 

Cp- specific heat 
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f - Coriolis parameter 
a - vertical coordinate (o = p/ir) 
a - vertical velocity (a = da/dt) 
q - Unir 

The continuity equation (1.4) may be refined by integrat- 
ing in the vertical and imposing boundary conditions on a; 
a(0) = 0(1) = 0. Introducing the notation G = V*Vq, we may 
rewr i te Eq . (1.4) as 

1-6) 1^ = -(D + G) . 

where the over bar denotes a vertical, average. 

The vertical velocity, a, may be obtained diagnostically 
by substituting Eq. (1.6) into Eq. (1.4) and integrating in 
the vertical to obtain 



-/ 



1.7) o = (0 + G)o - / (G + D) do 

o 

By defining a horizontally mean temperature. 



T = T (o) + T' (o. A, e, t) 



and using the following operator 



/V I 3X ^ 8Y 
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We can now expand the vector quantities in our 
equations in spherical, coordinates ( A , n ) to arrive 

I .8) = -a(A,B) 

1.9) = a(B,-A) -V^(r^E + (j) + RT%) 

I.IO) II- = -a(U0,ve) + 90 - 0 1^ 
dt 9a 

I . I I ) ||- = -(D + G) 

1.12) a 1^ = -RT 
9a 



where 

• 9V RT ^ 

A = (C+f)U ^ ~ ® 

r 

D - iTA.i\M A. ” RT' 9£n7T 

B - (c+f)V - o — - ^ _ 

r 

„ U 9£mr 3£n7r 

® = t:;t — * V -5j- 

r 

^ T~ 

2( I-U ) 

U = u cos n/r 
V = V cos n/r 
A = longitude 
n = latitude 
M = s i n n 



basic 

at 
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Equations (1.8) - (1.12) are the basic equations used in 
the model. The equations are represented spectrally in the 
horizontal and finite differenced in the vertical. The de- 
pendent variables are written in terms of a triangularly 
truncated series of spherical harmonics: 

X = z xJI y'J 



where the summation is a double sum over m,n for |m| _< M and 
|m| _< £ < M . The separation is such that the coefficients, 
X^, are functions of time and vertical coordinate and the 
spherical harmonics, Y^, are horizontal functions of space. 
The normalization and orthogonality properties of the Y^ 
allow the coefficients to be obtained as follows: 



X 



m 

£ 



<X, 




X y'J dydA 



The model used for this study has five vertical layers 
staggered as in Fig. I. The even levels carry a and the odd 
levels carry all the dependent variables plus the diagnostic 
variables, U,V. 

The finite difference scheme for the vertical advection 
terms is 



• _ I 

^ 3a 2(Aa, ) ‘-^^k+l 









This formulation conserves kinetic energy. 
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C7l = 0 





a^=0 






^,D.T,U.V 






07 


o^ 




4 ,d,t,u,v 






02 


^2 




^.D.T.U.V 

0-4 


-A 

03 




^^T,U.V 

05 


-N 
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VERTICAL SCHEME 



gure I. Staggerecj 5-Iayer sigma coordinate system. 
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Following Arakawa and Lamb (1976), the hydrostatic equa- 
tion (1.12) can be finite differenced as 

' • ^k+l ■ ^k ^ ^ ^'^k '''k+l ^k '''k^ 

where 

“k ' ' ®k ' ‘“k+/“k’'^^'^f - ' • 

A boundary condition is formed by integrating the hydrostatic 
equation from a = 0 to I which gives 




o o 



where is the terrain geopotential. This boundary con- 

dition is then finite differenced as 

1.14) Z (J)^(Aa^) - 4)^^^ = R Z T (Aa^) 
k 

Equations (1.13) and (1.14) are then combined to give the 
matrix equation 

1.15) 4) = [C]T + 4>3^^ 



which is the finite difference form of the hydrostatic equa- 
tion where the variables are now column vectors. 
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The pressure tendency equation (l.ll) can be written as 



at 



Z (G|^ + D|^) (Aa^) 
k 



or in matrix form as 



1.16) 



= -Cn] (G+D) 



Similarly equation (1.7) can be written in matrix form as 



1.17) 



a = [2] (G+D) 



The thermodynamic equation (I.IO) is differenced consis- 
tent with the method employed in the hydrostatic equation. 
The potential temperature is written in terms of temperature 
and the last term becomes 



ii 

3a 



2(Aa, 



f-. .^k >R/C^ 



k+l k+l - 

o 



k+l 



5 . R/C 
°k-I 



^^k - - 

^ (G. - G - D). 



For the purposes of the sem i - i mp I i c i t formulation de- 
scribed later in the paper, we separate the temperature into 
its horizontal average, T , and pe rtu rba t i on , T'. The 
matrix, Cy 3, is defined such that 












- T, 



T,. 

k ^ 



R/C. 



-) 



T 



k-l 



k+l 



2(Aa, ) 
k 






k-l 



2(Aa^) 



a, 



= CY]a 
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The vertically averaged divergence can be written D = 



CN]D, so that 





Cn]d . 



Substituting for a from equation (1.17) we obtain 
* * 

-CY]a - ^ cn:d = -CY]cz]G -(CYXz: + cnxd 

p p 

= -CS]G + CQ]D . 



Final ly the entire term can be written as 



, 18 ) 6 II = -CsMoDd- .-,^1 . .,^1 .. 



^ - T’) 



k+l 



rt; 



a. R/C RT. 

+ a, (t;-(^^^ — ) P t; , )] + 7 r-^(G. -G)- 7 ^ D 

k k ^ k-l C k C 

°k-i p p 



The nonlinear terms are computed using the transform 
method suggested by Orszag (1971). The longitudinal direc- 
tion is done using Fast Fourier Transform and the latitudinal 
direction is done with Gaussian Quadrature. The number of 
latitudes, N, and longitudes, M, satisfy: 

N ^ ( 3J+I )/2 and M ^ 3J+I . 

The number of points are picked to ensure non-a I iased results 
for quadratic terms. Certain terms containing a have triple 
products and are not computed al ias free. Experiments by 
the author and a previous study (Hoskins and Simmons, 1975) 
indicate that this source of error is negligible for the 
length of the integrations performed in this study. 
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111. SEMl-IMPLlClT SCHEME 



This section gives details of the sem i - i mp I i c i t time 
differencing scheme used in some of the experiments of this 
study. Using the matrices defined in the model description 
section, equations (I.9)-(1.I2) can be written 

2.1) 1^ + V^(4)' + RT%) = 

2.2) |I - [Q]D = F-^ 

2.3) 1^ + [N]D = -G 

2.4) ())' = CC]T 



where 



= (j) 



<P 



sfc • 



Fp and F^ represent the remaining terms in each equation 
which are not explicitly separated out. Following Robert 
et al. (1972), al I the terms on the left hand side of equa- 
tions (2.l)-(2.4) are time averaged. We define the following 
averaging operator and represent the local tendencies by a 
Leapfrog finite difference scheme: 






)t+At ^ 



( ) 



t-At 





( . ( )t-At 

2(At) 
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Equations (2.l)-(2.4) can then be written as 

2.5) 6^D = + RT^ 

2.6) 6^T = + CQ] d'*' 

2.7) 6^q = - G - CN] o"^ 

2.8) = CC] 

Eliminating (J) ’ ^ in equation (2.5) by substitution of Eq . 
(2.8), we arrive at a time averaged set of equations; 



2.9) 






+(At)(Fp - v^(CcJr^ 


2.10) 




jt-At 


+ (At) (F^ + CQ] d"^) 


2.11) 




t-At 

= q 


-(At) (G + Cn] D^) . 



Substituting equations (2.11) and (2.10) into (2.9) and 
solving for 0^ we obtain a Helmholtz equation for the time 
averaged divergence 

2.12) (DXAt)^[B]V^)D^ = + (At)CFj^-V^(CC]T''’"^’^ 

+ (At)Cc]F^+RT*q^~^''’-(At) RT^ G] 

where Cl] is the identify matrix and CB] = CC]CQ] ~ RT Cn]. 

For a spectral formulation, the equation is trivial to 

solve since the spherical harmonics are eigenfunctions of 

2 -t 

the Laplacian V . Knowing D , we then compute the prognostic 
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variables from equations (2.6) - (2.7). 

is computed from the definition of the time averaging 
operator : 






and the vort i c i ty 
forma I i sm and is 



Ca 1 cu 1 a t 


i ons 


with 


comp u te r 


1 1 me 


per 


possible 


with 


a t 



, never entered the sem i - i mp I i c i t 

computed from a standard Leapfrog scheme, 
the sem i - i mp I i c i t model required 5 % more 
time step. Stable integrations were 
ime step of 60 mins. 
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IV. QUAS I -GEOSTROPH I C DIVERGENCE 



Following the analysis of Simmons and Hoskins (1976), a 
diagnostic divergence is obtained for the case of quasi- 
geostrophic flow with the full variation of the Coriolis 



parameter over 


the sphere. The relevant equations become 


o.i; 


• V(?+f) - 5CD 



3.2) ^ - V^(({)' + RT*q) = 0 



3 3) — = - 

0.0) 


v’ . VT + a (^ I - |I ) 
V; C a 3a 

P 


where V , - 

•«> 

a 


rotational part of V 
sigma dot at odd levels 



= 7 • ( fV(V'^D) ) 



Making use 


of the hydrostatic relationship, equation 


(1.15), we ma y 


write equation (3.3) as 


d<p ' 
3t 


-CC](V^^ • VT) + y"" S 


where y = 


* * 
rci (— 1 _ 21 ) 

p 



2 

Taking V of this equation, aC of equation (3.1) and 
of eq . (3.2) and substituting we arrive at the fol lowing 



/N 

equation for d: 




3.4) y a = 


-v’^C (v^^-v(c+f) )+ ^d]+Cc](v^-vt)+rtV . 
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A second equation for a is obtained from the continuity 
equation in the form of eq. (1.17), which when explicitly 
wr i tten out. Is 



3.5) 



a = a(D + G) 




a 

(D + G) da 



where a is specified on even levels. 

/N 

The diagnostic a from equation (3.4) can be used to com- 
pute the divergence using eq. (3.5) provided we specify a 
boundary condition. For the purposes of initialization, this 
is most conveniently chosen to be 



9q _ 
at 



(D + G) 



0 



since this will tend to eliminate external gravity waves. 
Applying this constraint to equations (3,4) and (3.5) we ob- 
ta i n 

3.6) = -V‘^( (V^^-V(i;+f))+ S^O) t [C]V^*VT 

3.7) a = - I (D + G) da . 

o 

/N 

The method of solution will be to compute d using equa- 
tion (3.6) from an initial guess of the divergence. The even 
level d's are then found by interpolation and a new guess of 
the divergence is obtained by Inverting equation (3.7). The 
updated divergence can be substituted back into equation (3.6) 
and the procedure continues iteratively until sufficient 
accuracy is obtained. 



Equations (3.6) and (3.7) will now be put into spectral 
form. The equations contain three similar advective terms. 
It will be sufficient to show the spectral form of one of 
them; 



3.8) V^*V(C+f) = -J(i|^, ? +f) 



I fdllj dZ _ 9(5"*'f)\ 

"2 9l ■ 9l 



First, the terms are individually transformed to physical 
space by the summations: 



8p: . . 

m Jl ^imA 



_ 'ST' , 

3U " ^ 3y 
Z,m 

dZ . ^m _m ^imA 



Z,m 



. , m -m „ i mA 

^ = Pi *• 

J6,m 



d(Z+f) 

3y 



9P 

= + 2f! . 



£,m 



’2 9y 



The multiplications in eq . (3.8) are then computed at 

grid points and the result is transformed back to spectra 

space by computing the inner product defined as: 

+ 1 ,2tt 



3.9) 

L 4tt 




J(i|;, C+f)yJ^ dAdy = <J(i|), ?+f), Y^> . 



0 



The integral is computed using Gaussian quadrature in 
the y direction and Fast Fourier Transforms in the A direc- 
tion. The first term of equation (3.6) can then be 
eva I ua ted : 

SC(V^-V(?+f )) = V -(fV(V^J)) 

= fj -I- ^ ^ (V^J) 
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P-1 lon^ 

-V 

f J V . ' ’" ^* 1 t ■ M 1 f»D t 

f j T ^ t "^-.! ' 



-7a 



s. 



and expanded in spectral form using the previous computation 



of J 



< (Vi|;*V(C+f) ) y[^> 



+ 1 

= 



+ i 

Q I y(ZJ P^)pJ^dy-fi f (l-y^)(Z ^ 

-I 5- y. 



8P 



“ Ul+\) ay '■ L 



)P, dy , 



where the A integration has already been done. 
Making use of the recurrence relationship 



ypM = gM M ^ M M 

^^ 2 , ^£+1 ^£+1 ^£ £- 



M 






we can write the first integrand as 



E « j" p" = I e” p" + r-j" p” 

j ‘ ■^1-1 I j '^t+l ®«,+ l 



and using the orthogonality of the Legendre polynomials we 

obtain for the first integral 
-+l 



3.10) 



7 

-I 



,M ol^^oM ^ ,M M ^ ,M M 
y(Z P^)P^_ dy = J^_, . 



The second integral can be done similarly using the 
relationship 



3P^ 

, , 2, £ ^ M . M „M 

(l-y - (£ I) P^_| - ££j^^| 



26 



¥AL 







% 












I 






i 

I 



j| 



to obtain 



+1 ,M .pM 

'•">y ‘'-I' x^IaTnSiT'L 



M M 

5, dM , ,M ^L+l ,M 

= ^L+l — - \-l - 



Combining eqs (3.10) and (3.11) we obtain 



3.12) <4’(v^.v(c.f),y;;> = ncifi < ^ 



The second order operator X. D can now be formulated. 

,M. 

L' ' L 



Using the notation, <aCo, Y^> = , we can write by analogy 



with equation (3.12) 



pM _ „f-L+l , L ,M ^ 

Fl - a— D|_,| • 



Thus, 






Z 

I , m 



pM yM 

Z Z 



and 



vK - ^ - npL+l M ^ L M ^ 

X D, Y|^> if(SfD), Y^> nC E|_ F|^_| + '^l+I^ 



M M 

Substituting in for F, , and F, .,, we obtain the result 
^ L- I Lt I ' 



3.13) < 



-<^D,Y^> = 



^2r,(Lt| )L 

^ ^ L(L-I ) L L-l °L-2 

^ ,L^-I , M,2 , L^+ 2L , M ,2,^M 

+ ( (e^) + )D^ 



+ ( 



+ L 



( L+ 1 ) ■ 
M M 



(L+KL+2) "L+l L+2' "L+2 



X rS 
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Finally, eqs (3.6) and (3.7) may now be fully expanded 



as 



3.14) + ^^.^^(<<j((V^-VC+f),Yj^> + <rf^D,vJ^>)] 



T ICS - /* 4. r^^ 

3.15) ~ ~ J 



where 



Rl = He] <v^ • VT, Y^> 
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V. DESIGN OF THE EXPERIMENTS 



The initial conditions used in the experiments of this 

thesis were developed from a baroci i n i ca I ly unstable mean 

flow and a small perturbation in zonal wave no. 6 (m=6). 

The mass field of the mean flow was constrained to be in 

geostrophic thermal balance with the mean wind. The vertica 

prof i le of the mean wind was I inear varying from 4.5 m/sec 

at O’ = I to 49.5 m/sec at ct = 0. The latitudinal structure 
2 

varied as sin (2r\) placing a jet max at n = 45°. The verti- 
cal temperature profile approximated the U. S. Standard 
Atmosp here . 

The perturbation in wave no. 6 corresponded to a maximum 
v-component of .1 m/sec. The model was integrated for 192 
hrs allowing the baroclinic wave to grow to finite amplitude 
Fig. 2 shows a time series of the vorticity during this 
period. The growth is strikingly linear until 192 hrs where 
there is evidence of the mean flow (m = 0) changing due to 
feedback from the eddies. The I inear growth corresponds to 
a doubling time of roughly 26 hrs for m = 6. The fields at 
192 hrs provided the balanced fields for the initialization 
exper iments . 

Three types of experiments and one control case were 
integrated for 72 hrs of model time. 
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Figure 2. Vorticity time series during the first 192 hrs. 

Vertical axis is logarithmic. Units are secs"*. 
Solid curve is the magnitude of wave no. 0 at 
latitude 60°N. Dotted curve is m=6, short 
dashes are m=2 and long dashes are m=l8. 
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The control case consisted of initializing with the 
balanced fields and integrating. The experiments consisted 
of the following: 



Experiment A: 
Experiment B : 
Experiment C: 



Initialize with the balanced fields but 
no divergence 

Initialize with the balanced fields plus 
the qua s i -geostrop h i c divergence 

Identical to Experiment A but with a semi- 
imp I i c i t scheme 



All experiments and the control were first integrated 
with a spectral truncation of wave no. 21. These will be 
referred to as low resolution experiments. High resolution 
experiments were integrated for the control and Experiments 
A and B. These cases were truncated at wave no. 42. Time 
steps of 12 mins for the low resolution and 6 mins for the 
high were used. Experiment C was integrated with time steps 
of 12 and 60 mins. Al I experiments had a weak Robert time 
filter of 0.1 and no dissipation term was used. 
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V I . RESULTS 



A. CONTROL CASE 

Figures 3-6 show the initial fields at level 5 for the 
low resolution experiments. Figures 7-10 are the correspond- 
ing high resolution fields. They are cha racter i zed by a 
developing frontal zone. The advantage of the high resolu- 
tion is apparent in the temperature gradients across the 
cold front and warm front and in the surface pressure. The 
low in the high resolution surface pressure field is 5 mb 
deeper. The negative and positive areas of the divergence 
and vorticity fields are approximately symmetrical in the 
low resolution fields and quite asymmetric in the high reso- 
lution fields. The region of maximum positive vorticity 
corresponds very well with the maximum temperature gradient. 

The mean fields consist of an indirect cell induced by 
the growing baroclinic eddies. The indirect cell develops 
in order to decrease the vertical wind shear which must remain 
in approximate thermal balance with the decreasing pole to 
equator temperature gradient. The effect is more pronounced 
in the low resolution fields. The high resolution mean 
fields contain a weak indirect circulation and a strong 
direct cell. 

The dynamic situation during the period of the control 
case is one of a rapid cascade of energy to sma I I er scales. 
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Figure 3. Initial low resolution surface pressure. Contour interval 
is 5 mb . 
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Figure 5. Initial low resolution vorticity for o = .9. Contour 
Interval Is 5xI0"°sec-l. Mean field is scaled by 1 0^ 
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Figure 6. Initial low resolution divergence for a = ,9, Contou 
interval is IxIO^^sec"!. Mean field Is scaled bv 1 0^ 
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Initial high resolution vorticity for a = o9» Contou 
interval is 5xlO“^sec“l. Mean field is scaled by 10^ 



u 




-1 



it 



- i 



•J 



4 














H 



Ml"! 



"’lira 




4 



PERTURBRTION FIELDS 




LU 

a 



L. 

3 • 

oo 
-I- o 
c — 
o 
o 



O "O 
CJv 0 

o 

(U 

II u 
D 

CO 

u — 

o 

TU 

0 0 

u — 

C H- 
0 

O) c 
u (u 
0 0 
> 2: 




DIVERGENCE IN 1/SECS 



c — 

O I 

— o 

-I- 0 
3 lO 
^kO 
O I 
(O O 
0 _ 
^ X 



U) 10 
JZ 



— (U 
(U > 

— u 
+- 0 

— -h- 
C C 



o 



0 

u 

3 

O) 



40 



The horizontal scale of the frontal zone decreases. Figures 
11-13 show this effect in the high resolution temperature 
fields for times out to 36 hr. Beyond this point, the 
energy builds up at the spectral limit of the model since 
there is no dissipation and the 72 hr fields are non-meteoro- 
logical (Fig. 14 and 15). Because the energy blockage occurs 
at smaller scales, the high resolution case (Fig. 14) is, 
paradoxically, noisier and in greater error. 

B. GRAVITY WAVES 

In anticipation of later sections, we discuss the 
character i st i cs of inertial gravity waves for the model used 
in this study. The divergence is the most sensitive indica- 
tor of gravity wave motion. As discussed in Hoskins and 
Simmons (1975), equation (2.12) is the sem i - i mp I i c i t analogue 
of the gravity wave equation for a muiti~layer model. The 
eigenvalues of the matrix Cb 3 (eq. 2.12) give the permissible 
e i gen ve I oc i t i es of the gravity wave modes. The corresponding 
eigenvectors give the vertical structure. Figure 16 depicts 
the eigenmodes for the experiments in this study. Table I 
gives the periods of the gravity wave modes for zonal wave 
no . 6 in h rs . 



TABLE I 



Exte rna 1 


First 


Second 


Th i rd 


Fou rt h 


3,3 1 


1 1 . 1 


28.0 


61.7 


1 54 
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Figure 15. 72 hr low resolution temperature for a = .9o Contour 
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internal divergence modes* Velocities are in m/sec« 



The total divergence may be thought of as consisting of 
a sum of the individual eigen modes: 

4.1) D = S D|^ 

k 

where the are the eigenvectors. The weighting coefficients, 
a^, may be obtained from the simultaneous equations formed 
by the dot products of the individual eigenvectors with equa- 
tion (4.1). A time series of the coefficients will provide 
useful diagnostic information in later sections. Figures 17 
and 18 are a time series of the external and first internal 
modes for the low resolution control case at a point on the 
60° latitude circle. It is apparent that the coefficients 
contain both the meteorological modes and any excited gravity 
modes. Figure 17 contains no high frequency oscillation and 
is purely meteorological. Figure 18 contains a meteorologi- 
cal mode with a period of roughly 3 days plus a sma I ler 
amplitude wave whose period of approximately 1/2 day corre- 
sponds well with the calculated value for the first internal 
mode from Table I. The remaining internal gravity modes 
have periods approaching meteorological values and are diffi- 
cult to distinguish. 

C. EXPERIMENTS A AND B 
I . Low Resolution 

Neglecting the divergence in the initialized fields 
(Experiment A) generated computational gravity waves. Compar- 
ing the pressure tendencies of the control case (Fig. 19c) 
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Figure 17. Divergence time series of the external 
mode for the low resolution control 
case. Vertical axis units are sec~^. 
Time is in hrs. 
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Figure 18. 



Time series of the first internal mode 
- otherwise as in Fig. 17. 
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with those of Experiment A in Fig. 19a, it is apparent that 
high frequency oscillations of roughly 3 hr period are 
excited. The initial magnitudes of the oscillations are 
approximately 1-2 mb/hr. The oscillations dampen out with 
time as the geostrophic adjustment takes place although the 
dispersal of the gravity waves during the adjustment process 
is hampered by the longitudinal periodicity and the symmetry 
of the fields about the equator. 

Figure 19b shows the pressure tendencies with the 
q ua s i -geos t rop h i c divergence inserted into the initializa- 
tion fields (Experiment B). The highest frequency oscilla- 
tions have been reduced ove r . Exp e r i men t A by a factor of 
a p p rox i ma te I y 2 . 

Comparison of the divergence modes in Figs 20, 21, 22 
allows us to be more quantitative. The figures have the 
control case subtracted out and contain the gravity wave 
noise generated by the approximations of each experiment. 

The noise generated in the internal modes by the non -d i ve r gen t 
initialization is an order of magnitude larger than the ex- 
ternal mode and is generally 25^ of the total divergence. 

The figures show that the external mode of Experiment B is 
about 68^ of Experiment A. The first internal mode shows 
more noise in Experiment B initial ly and then decreasing to 
50^ and the second internal mode shows a decrease of about 
50^ over Experiment A. In general the effect of the quasi- 
geostrophic divergence is to reduce the gravity wave noise. 
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Figure 20. Gravity wave noise in the external mode for low resolution 
Experiment A (20a) and Experiment B (20b). Vertical scale 
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Figure 21. Gravity wave noise in the first Internal mode - otherwise 
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High Resolution 



Figure 23 shows the quas i -geostroph i c divergence at 
level 5 which was used for Experiment B. A qualitative com- 
parison with the actual model divergence (Fig. 10) shows the 
phases to be in agreement but that the q ua s i -geost rop h i c 
approximation is an overestimate for the perturbation fields. 
The zonal mean fields do not compare well. The quasi- 
geostrophic divergence shows a strong indirect circulation 
which is displaced poleward relative to the model divergence. 
There is no evidence of a direct cell. 

Figures 24-26 depict the noise generated in the 
divergence modes for Experiments A and B. General ly they 
show that the gravity wave noise of Experiment B is 66 %- 12 % 
of that in Experiment A. The effect of the q ua s i -geo st rop h i c 
divergence is to reduce the noise. However, this reduction 
is less than in the low resolution experiments. 

In view of the fact that the high resolution experi- 
ments were less successful than the low resolution and that 
qu a s i -geost rop h i c formulations are valid only for small 
values of the Rossby number, it is interesting to compare the 
initial model divergence of the control case with the quasi- 
geostrophic divergence in a spectral manner. Referring to 
Fig. 27, it is seen that there is good agreement between the 
model and qua s i -geost rop h i c divergences for zonal waves 
cor respond i ng to m = 6 , 12, 18. The overestimation by the 

quas i -geost roph i c divergence noted earlier is primarily in 
wave no. 6. The mean fields (m = 0) and waves correspond i ng 
to m = 24,30,36,42 do not compare well. 
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Figure 24. Gravity wave noise in the external mode for high resolution 
Experiment A (24a) and Experiment B (24b). Vertical scale 
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Figure 25. Gravity wave noise in the first interval mode - otherwise 
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Figure 26. Gravity wave noise in the second Interval - otherwise 
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Figure 27. Initial divergence coefficients for the control case 
(dotted line) and Experiment B (dashed line). Verti 
axis is latitudinal wave number (lat. index of the 
spherical harmonics). Horizontal scale is sec” ^ . 



One motive for initializing with the quas i -geost ro- 
phic divergence is that, presumably, the "true" divergence 
would be approximated at the start and that the model would 
thus develop this divergence field earlier than a non- 
divergent initialization. Figures 28 and 29 address this 
point. They show that at 6 hours into the integrations, zonal 
wave numbers 6, 12, L8 do, in fact, approach the control 

divergence slightly faster. There is almost no difference 
for the higher wave numbers and the mean field shows poor 
adjustment in either case. This might be expected for the 
mean field since the geostrophic adjustment process is very 
inefficient at this scale in a bounded region. On the other 
hand, the higher wave numbers adjust very rapidly and after 
6 hrs, they are in good agreement with the control case. 

In light of the results showing good agreement in 
the initial quas i -geost roph i c divergence with the control 
divergence for zonal wave number 6, 12, 18 and bad agreement 

for other scales, it is natural to modify Experiment B such 
that only scales corresponding to m = 6, 12, 18 are used 

in the quas i -geostrop h i c divergence. Figure 30 shows the 
time series of the external divergence mode for this experi- 
ment. Figure 31 shows the spectral adjustment of the diver- 
gence at 6 hours. The corresponding results for the unmodi- 
fied Experiment B are Figs. 24 and 29. The results are 
almost identical to the unmodified quas i -geostroph i c diver- 
gence results. 
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Figure 28. 6 hr divergence coefficients for the control case (dotted 

line) and Experiment A (dashed line) - otherwise as in 
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Figure 29. 6 hr divergence coefficients for the control case (dotted 

I ine) and Experiment B (dashed I ine) - otherwise as in 
Fig. 27. 




TIHE 



Figure 30. Time series of the gravity wave noise 

in themodified Experiment B. Vertical 
scale is sec“l. Time is in hrs. 
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Figure 31. 6 hr divergence coefficients for the control case (dotted 

line) and modified Experiment B (dashed line). Vertical 
axis latitudinal wave no. Horizontal scale is sec“*. 



D. EXPERIMENT C 



As discussed in Hoskins and Simmons (1975), a linear sta- 
bility analysis of the sem i - i mp I i c i t scheme shows that the 
time step is limited only by the Rossby waves. The effect 
of the scheme is to slow down the faster moving gravity wave 
modes so that the integrations do not amplify in time. 

Slower moving atmospheric motions are less affected. Table I 
presents sample calculations of the distortion of the periods 
of the gravity waves present in these experiments for two 
time steps. These may be compared with Table I. The general 
behavior is that the periods increase with increasing time 
steps. 



TABLE I I 



At 


Externa 1 


First 


Second 


Th i rd 


Fourth 


12 mins 


3.31 


11.14 


28.0 


6 1.7 


154.0 


60 mins 


5.7 


12.2 


28.5 


6 1.9 


154.1 



Two points will be investigated in this experiment. The 
first is that the sem i - i mp I i.c i t scheme should tend to smooth 
out high frequency oscillations, particularly with a large 
time step. The second is that the sem i - i mp I i c i t method 
might possibly hamper the geostrophic adjustment process 
since it slows down the gravity waves which are the mechanism 
for the adjustment. 
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Figure 32 shows the pressure tendencies using the semi- 
implicit method with time steps of 12 and 60 mins. Compar- 



ing this with the corresponding results for the explicit 
integrations (Fig. 19a), it is evident that the effect of 
the sem i - i mp I i c i t method is to smooth out the 3 hr oscilla- 
tions. Figure 33 shows the relative effect of a small vs. 
large time step in the semi-impi iclt method on the external 
gravity waves. The larger period of the 60 min time step is 
evident. The larger time step initially amplifies the com- 
putational noise, but after 24 hrs it is considerably damped, 
The smoothing in the pressure tendencies is postulated to 
come about from the decreased velocity of the external 
gravity mode assuming an advective time scale. 

The second point of this experiment is shown in Figs. 34 
and 35. Comparison of the divergence coefficients for time 
steps of 12 and 60 mins shows no difference at 3 hrs. The 
adjustment of the divergence to the control case is apparent- 
ly the same for either time step. 



I 



68 




‘ i 



■M 



u 



i I#' 

^'•1 




t 

fH 



: I * 

-O .fHi 




u 

• \ 

e 




•Si 






iJE-06 



I 





CM 




69 



Figure 32. aq/3+ for Experiment C with a time step of 12 mins (32a) 

and 60 mins (32b). Verticai scale is sec"l 
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Figure 33. Gravity wave noise in the externai mode for Experiment 
with time steps of 12 mins (33a) and 60 mins (33b). 
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Figure 34 . 



3 hr divergence coefficients for the 
resolution control case (dotted line) 
Experiment C with a 12 min time step 
(dashed line). Vertical axis is lati 
wave no. Horizontal scale is sec"l. 
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hr divergence coefficients with a 60 
n time step - otherwise as in Fig. 34. 
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VII. SUMMARY AND CONCLUSIONS 



The results of this study are based upon the divergence 
developed by a baroclinic wave on a uniform earth. Initiali- 
zation without the divergent wind component generated spurious 
inertial gravity waves. The internal state of the initial 
fields had nearly zero mean divergence. Thus, the external 
mode of the generated gravity waves was an order of magni- 
tude less than the total divergence and was not significant. 
However, the internal gravity wave modes generated were signi- 
ficant. An important source for the external mode, namely 
topography, was neglected. 

Initializing with the q ua s i -geost rop h i c divergence 
generally decreased the spurious gravity wave noise and aided 
the ensuing integrations to develop a divergence closer to 
the control case, but both effects were small. Furthermore, 
the advantage of the q ua s i -g eos t rop h i c divergence decreased 
when smaller scales were included in the integrations. In- 
serting the qua s i -geostrop h i c divergence into only synoptic 
scales gave no advantage. However, the effect of this experi- 
ment may have been reduced by the fact that there were no 
planetary scale waves. 

Using a sem i - i mp I i c i t scheme with the same time step as 
an explicit integration had no effect on the divergence modes 
but did smooth out oscillations in the pressure tendencies. 
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A larger time step initially increased the amplitude and 
period of the fast external gravity waves and had less effect 
on slower moving motions while still stabilizing pressure 
tendencies. 

In general, the effect of a divergent initialization for 
a global primitive equation model is small. This same con- 
clusion has been reached by two similar studies - Houghton 
et al (1971) and Dey et al (1975). This study, however, has 
neglected any physical processes such as precipitation as 
wel I as orography and the effect of an initial divergence in 
these cases may well be more important. 
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